% clear all
% %clc
% 
% load('MC_N10K_m3_b4_rand_1_fullinfo.mat')

load('grid_100.mat')
N_temp = 10000;

% [data,data_trans] = generate_data_rand(N_temp,J,mu_x,2,mu_z,2,alpha,beta,fraction);

slope = zeros(nMC,1);
ratio = zeros(nMC,1);
ratio_means = zeros(nMC,1);
ratio_medians = zeros(nMC,1);
ratio_trimmed = zeros(nMC,1);
ratio_trimmed_mean = zeros(nMC,1);

for ss=1:nMC
    
    ratio_temp = zeros(N_temp,1);
    ds_dz = zeros(N_temp,1);
    ds_dx = zeros(N_temp,1);
    
    parfor i=1:N_temp
               
        argument = normcdf([data.x(i,:) data.z(i,:)]);
        
        ds_dz(i) = fun_dsigma(theta_NPD(:,ss),argument,m_regr,combos_all_unique,J+1)*normpdf(data.z(i,1));
                
        ds_dx(i) = fun_dsigma(theta_NPD(:,ss),argument,m_regr,combos_all_unique,1)*normpdf(data.x(i,1));
        
        ratio_temp(i) = ds_dz(i)/ds_dx(i);
        %i
    end
    
    ratio(ss) = median(ratio_temp);
    
    ratio_means(ss) = mean(ds_dz)/mean(ds_dx);
    ratio_medians(ss) = median(ds_dz)/median(ds_dx);
    
    [num_temp_ord,ind_num] = sort(ds_dz);
    [den_temp_ord,ind_den] = sort(ds_dx);
    ratio_trimmed_mean(ss) = mean(num_temp_ord(N_temp*0.25:N_temp*0.75))/mean(den_temp_ord(N_temp*0.25:N_temp*0.75));
    ratio_trimmed_mean(ss)
    
    [ratio_temp_ord,ind] = sort(ratio_temp);
    ratio_trimmed(ss) = mean(ratio_temp_ord(N_temp*0.10:N_temp*0.90));
    
    regr_coeff = regress(ds_dz,ds_dx); % regr_coeff = regress(d2s_dzdz,[ones(length(d2s_dzdz),1) d2s_dzdx])
    
    slope(ss) = regr_coeff(1); %regr_coeff(2);
    
    ss
    
end
